Introducción a R para la Economía y el Análisis de Políticas Públicas

Desde La Mancha con amor

Autor/a

Roberto Martínez Lacoba (@robertomlacoba en ciertas redes)

Fecha de publicación

24 de enero de 2025

0. Manuales interesantes

En este apartado os dejo información sobre algunos manuales y recursos que han sido utilizados para el desarrollo de este documento, así como otros que me han acompañado durante los últimos años.

Introducción a la programación en R

Libros de econometría y estadística aplicada

Gráficos y visualización

Epidemiología para salud pública

Microeconometría y análisis económico

Bioestadística y ciencia de datos

Otros

  1. R for Data Science Este es un clásico joven que han hecho los creadores del tidyverse y también tienen cositas en ggplot2 - gente fiera, vaya.
  2. Advanced R
  3. Hands-On Programming with R

1. Recordatorio de instalación de R, RStudio y Quarto

  1. Instalación de R y RStudio: https://posit.co/download/rstudio-desktop/
  2. Instalar Rtools y devtools: Para descargar paquetes de github
  1. Instalar Quarto: https://quarto.org/docs/get-started/

Paquetes que (quizá) utilizaremos

paquetes <- c("GGally", "HistData", "RColorBrewer", "caret", "CDR", "classInt", "corrplot","data.table", "DiagrammeR", "devtools", "DT", "eurostat", "explore", "factoextra", "gapminder", "ggalluvial", "gganimate", "ggmosaic", "ggplot2", "ggpubr", "ggrepel", "ggridges", "ggstatsplot", "gplots", "gstat", "ggthemes", "igraph", "knitr", "leaflet", "leaflet.extras", "leaflet.providers", "lmtest", "mapSpain", "openxlsx", "pacman", "pROC", "plotly", "psych", "quarto", "raster", "readxl", "remotes", "rlang", "rvest", "reshape2", "scales", "sf", "sjPlot","skimr", "summarytools", "survival", "tidyverse", "titanic", "tm", "visdat", "viridis", "wordcloud", "quantmod","fBasics","curl","dygraphs","tidyquant", "timetk","highcharter", "wooldridge")

instalarpaquetes <- function(x){ # Creamos nuestra primera función :)
  for (i in x){ # Para cada ítem en el objeto x
    if (!require(i, character.only = TRUE)){
      install.packages(i, dependencies = TRUE)
      library(i, character.only = TRUE)
    }
  }
}

instalarpaquetes(paquetes)

Antes de avanzar, recuerda es relevante fijar el directorio de trabajo. Si estás en un proyecto de R o de Quarto el directorio se ubica en ese mismo proyecto. Si simplemente trabajas con carpetas en Windows (no tiene por qué estar mal hacerlo así), recuerda fijarlo en la carpeta correspondiente.Se puede predefinir en las Global Options de RStudio.

getwd() # Muestra el directorio de trabajo
[1] "C:/Users/Roberto.MLacoba/OneDrive - Universidad de Castilla-La Mancha/UCLM/Clases universitarias/Lo básico en R/Intro R Eco y APP"
#setwd() # Fija el directorio de trabajo en la ruta asignada, por ejemplo: C:/Users/Usuario/Desktop

2. Operadores elementales e iniciación a la programación

2.1 Operadores elementales

2.1.1 Operadores aritméticos

  1. Suma: +
  2. Resta: -
  3. Multiplicación: *
  4. División: /
  5. Exponenciación: ^
  6. Módulo: %% (resto de la división)
  7. División entera: %/% (cociente de la división)
  8. Raíz cuadrada: sqrt()
  9. Logaritmo natural: log()
  10. Logaritmo en base 10: log10()
  11. Exponencial: exp() e^x
  12. Valor absoluto: abs()
  13. Factorial: factorial()
  14. Producto de matrices: %*%

# Ejemplos de operadores aritméticos

# Ejemplos de operadores aritméticos
x <- 5
y <- 3
z <- 2

x+y # Suma
[1] 8
x-y # Resta
[1] 2
sqrt(y)
[1] 1.732051
exp(1)
[1] 2.718282
log(10)
[1] 2.302585
log10(100)
[1] 2
log(exp(1))
[1] 1
a<- as.matrix(c(1,2,3)) #3x1
b<- as.matrix(c(4,5,6)) #3x1
bt<-t(as.matrix(c(4,5,6))) #1x3)
a%*%b # No funciona
Error in a %*% b: argumentos no compatibles
a%*%bt # 3x1 * 1x3 = 3x3
     [,1] [,2] [,3]
[1,]    4    5    6
[2,]    8   10   12
[3,]   12   15   18
bt%*%a # 1x3 * 3x1 = 1x1
     [,1]
[1,]   32

2.1.2 Operadores de comparación

  1. Mayor que: >
  2. Menor que: <
  3. Mayor o igual que: >=
  4. Menor o igual que: <=
  5. Igual que: ==
  6. Distinto de: !=

2.1.3 Operadores lógicos

  1. Y: &
  2. O: |
  3. No: !
  4. Y exclusivo: xor()
  5. Todo verdadero: all()
  6. Al menos uno verdadero: any()

2.1.4 Operadores de asignación

  1. Asignación simple: <- o = o ->
  2. El pipe %>% (de magrittr) o también |> (de base R): se puede usar para encadenar funciones y se interpreta como un “entonces”. Traducción literal: tubería

Vemos unos ejemplos con 2.1.2, 2.1.3 y 2.1.4 utilizando R base.

# Usamos la clásica base de datos iris
a<- mean(iris$Sepal.Length[iris$Sepal.Length>=5]) # Media de la longitud del sépalo si es mayor o igual que 5
a
[1] 6.041406
b=mean(iris$Sepal.Length[iris$Sepal.Length>=5 & iris$Species=="setosa"]) # Media de la longitud del sépalo si es mayor o igual que 5 y la especie es setosa
b
[1] 5.23
mean(iris$Sepal.Length[iris$Sepal.Length>=5 | iris$Petal.Length>=1.5]) ->c # Media de la longitud del sépalo si es mayor o igual que 5 o la longitud del pétalo es mayor o igual que 1.5
c
[1] 5.960584

Como se observa, R ofrece una versatilidad razonable ya que permite subdividir de forma sencilla cualquier muestra/información. Imagina cómo harías esto en Excel… Es más, vamos a hacerlo y comparamos (y de paso aprendemos a exportar documentos en formato xlsx).

# Exportamos a Excel
library(openxlsx)
write.xlsx(iris, "data/iris.xlsx")

2.1.5 Operadores de concatenación

  1. Concatenación de caracteres: paste() y paste0() 
a <- paste("Hola", "Roger:", sep=" ")
b <- paste("tot", "bé?", sep=" ")
paste(a, b, sep=" ")
[1] "Hola Roger: tot bé?"
  1. Concatenación de vectores: c()
  2. Concatenación de listas: list()
  3. Concatenación de matrices y vectores: rbind() y cbind()

2.2 Iniciación a la programación

Programar es un arte, un arte divertido, pero nadie nace aprendido y el ensayo y el error es clave (esto también aplica al manejo de R). En este apartado se recogen algunos elementos fundamentales para que sirva de introducción: tómalo como si fuese un juego. Puedes intentar reproducir el código a mano (así practicas escribiendo a mano) o cópialo, úsalo y altéralo para ver qué pasa. Lo que sí es un aspecto relevante de la programación (y no solo de la programación, también en otras disciplinas) es la posibilidad de imaginar. Cuando queramos hacer algún tipo de alteración de una base de datos debemos pensar qué es lo que realmente queremos hacer: ¿queremos obtener algo que responde de condiciones?, ¿queremos reproducir y repetir una serie de valores siguiendo un determinado patrón?, ¿queremos que ocurra algo si se cumple una condición? Pensemos, en este sentido, como con el lenguaje natural. Decía un filósofo, no recuerdo el nombre, que el mundo de las personas es tan amplio como las palabras que conocen. En estas disciplinas pensemos de forma similar: podremos hacer tanto como podamos ser capaces de imágenes; seguir formándonos y aprendiendo es fundamental. El diablo sabe más por viejo que por diablo.

2.2.1 Estructuras de control

  1. Condicionales: if, else, else if
  2. Bucles: for, while, repeat
  3. Interrupción de bucles: break, next
  4. Funciones: function()
  5. Estructuras de control: switch(), tryCatch(). Esto sirve para evaluar una opción (switch) o para capturar errores (tryCatch).

Profundizaremos en los condicionales, en los bucles y en su interrupción.

2.2.2 Ejemplos de condicionales, bucles y de interrupción de bucles

Partimos de unos ejemplos sencillos sobre la propia consola/script de R con el fin de familiarizarnos con R y su entorno:

# Ejemplo de condicional
 x <- 5
 if (x > 5){
   print("x es mayor que 5")
 } else {
   print("x es menor o igual que 5")
 }

 # Ejemplo de bucle
 for (i in 1:5){
   print(i)
 }
 # Ejemplo de bucle creando un vector:
 vector <- c() # Definimos un vector vacío
 for (i in rep(1:10,times=2, each=5)){
   vector <- c(vector, i) # Lo que hace es iterar y añadir el valor de i al vector
 }
 print(vector)
 
 vector <- c()
 for (i in seq(1,10,2)){
   vector <- c(vector, i)
 }
 print(vector)
 
 # Ahora (muchísimo) más eficiente:
 vector <- rep(1:10, times=2, each=5)
 
 # Nota: pensar siempre cuál es la mejor forma de llegar al resultado
 
 # Ejemplo de bucle con interrupción
 for (i in 1:10){
   if (i == 5){
     break
   }
   print(i)
 }
 # Ejemplo con repeat:
 i <- 1 # Número de partida
 repeat{ # Repite
   print(i) # Imprime i
   i <- i + 1 # Suma 1 a i
   if (i > 100){ # Si i es mayor que 100 
     break # Para
   }
 }
 
 # Ejemplo con while:
 i <- 1
 while (i <= 10){ # Mientras i sea menor o igual que 10
   print(i) # Imprime i
   i <- i + 1 # Suma 1 a i 
 }
 
 # Ejemplo combinando bucles y condicionales
 for (i in 1:10){
   if (i %% 2 == 0){ # Si el resto de dividir i entre 2 es 0
     print(paste(i, "es par")) # Imprime i y que es par
   } else { # Si no
     print(paste(i, "es impar")) # Imprime i y que es impar
   }
 }
 
 # Ejemplo con next:
 for (i in 1:10){
   if (i %% 2 == 1){ # Si el resto de dividir i entre 2 es 1
     next # Siguiente
   }
   print(i) # Imprime i
 }

2.2.3 El uso de condicionales y bucles sobre bases de datos

Para realizar esta práctica vamos a utilizar los datos que usé para una investigación pasada que está publicada en este enlace. Preparamos el chunk de R para cargar los datos (recuerda que para cargar los datos en R hay que seleccionar una ruta determinada):

# Cargamos los datos
datos <- openxlsx::read.xlsx("data/Dataset.xlsx", sheet=2) # Otra forma de "invocar" funciones de paquetes sin pasar por library. Recuerda **adaptar tu ruta**
head(datos) # Para ver las primeras filas. Se puede configurar como head(x, 3), etc.
  ID age gender hiseilevel location cook degree    dairy olivesnuts herbsspices
1  1  20      1          1        1    0      0 3.421918  0.7123288  0.27945205
2  2  45      1          2        1    1      0 6.643836  2.0000000  2.00000000
3  3  20      0          3        1    0      1 3.421918  0.4273973  0.00000000
4  4  18      1          3        1    0      0 1.498630  0.1315068  0.06575342
5  5  18      0          1        2    1      0 1.263014  0.0000000  1.21369863
6  6  18      1          2        1    0      0 9.627397  0.2794521  0.27945205
     fruit    veget   oliveoil breadpasta   potato   redmeat     sweets
1 4.115068 1.397260 0.49863014  2.0684932 1.495890  5.753425  6.3671233
2 4.279452 7.852055 2.50136986  5.1424658 0.000000 11.832877  1.9561644
3 1.463014 1.827397 2.50136986  5.2082192 0.460274  8.860274 18.8712329
4 1.895890 2.608219 0.49863014  1.6904110 0.460274 12.312329  8.8986301
5 3.610959 2.936986 0.06575342  0.9863014 0.460274  5.753425  0.9205479
6 5.117808 2.043836 1.00000000  2.0684932 1.495890 20.328767  4.2575342
  whitemeat     fish      egg   legumes alcoholicUA  fastfood precooked
1 1.4958904 5.408219 1.495890 3.4520548       0.210 0.3287671 1.7643836
2 4.9863014 3.452055 4.986301 5.4082192       0.000 0.2630137 2.5671233
3 1.9561644 3.912329 1.495890 2.8767123       0.000 0.9205479 1.2547945
4 1.9561644 4.717808 3.490411 2.3013699       0.000 0.4602740 1.5260274
5 0.9205479 9.819178 0.460274 0.9205479       3.382 0.7068493 1.0246575
6 6.4821918 6.789041 3.490411 2.4164384       0.000 0.7616438 0.9041096
str(datos)
'data.frame':   593 obs. of  24 variables:
 $ ID         : num  1 2 3 4 5 6 7 8 9 10 ...
 $ age        : num  20 45 20 18 18 18 18 18 18 18 ...
 $ gender     : num  1 1 0 1 0 1 0 0 1 0 ...
 $ hiseilevel : num  1 2 3 3 1 2 3 1 2 2 ...
 $ location   : num  1 1 1 1 2 1 2 3 1 1 ...
 $ cook       : num  0 1 0 0 1 0 0 1 0 0 ...
 $ degree     : num  0 0 1 0 0 0 0 0 0 0 ...
 $ dairy      : num  3.42 6.64 3.42 1.5 1.26 ...
 $ olivesnuts : num  0.712 2 0.427 0.132 0 ...
 $ herbsspices: num  0.2795 2 0 0.0658 1.2137 ...
 $ fruit      : num  4.12 4.28 1.46 1.9 3.61 ...
 $ veget      : num  1.4 7.85 1.83 2.61 2.94 ...
 $ oliveoil   : num  0.4986 2.5014 2.5014 0.4986 0.0658 ...
 $ breadpasta : num  2.068 5.142 5.208 1.69 0.986 ...
 $ potato     : num  1.5 0 0.46 0.46 0.46 ...
 $ redmeat    : num  5.75 11.83 8.86 12.31 5.75 ...
 $ sweets     : num  6.367 1.956 18.871 8.899 0.921 ...
 $ whitemeat  : num  1.496 4.986 1.956 1.956 0.921 ...
 $ fish       : num  5.41 3.45 3.91 4.72 9.82 ...
 $ egg        : num  1.5 4.99 1.5 3.49 0.46 ...
 $ legumes    : num  3.452 5.408 2.877 2.301 0.921 ...
 $ alcoholicUA: num  0.21 0 0 0 3.38 ...
 $ fastfood   : num  0.329 0.263 0.921 0.46 0.707 ...
 $ precooked  : num  1.76 2.57 1.25 1.53 1.02 ...
# Arreglamos los datos de las primeras columnas para que sean factores y no numéricos
datos$gender <- as.factor(datos$gender)
datos$hiseilevel <- as.factor(datos$hiseilevel)
datos$location <- as.factor(datos$location)
datos$cook <- as.factor(datos$cook)
datos$degree <- as.factor(datos$degree)

Quizá se puede hacer de otra forma:

datos <- read.xlsx("data/Dataset.xlsx", sheet=2)
# Selección de columnas
columnas_a_factor <- c("gender", "hiseilevel", "location", "cook", "degree")

# Aplicando for
for (col in columnas_a_factor) {
  datos[[col]] <- as.factor(datos[[col]])
}

O de otra forma (de las más utilizadas actualmente):

datos <- read.xlsx("data/Dataset.xlsx", sheet=2)
library(dplyr)

datos2 <- datos %>% # Nuestro primer pipe :) -> Pongo datos2 por tener otra base distinta
  mutate(across(c(gender, hiseilevel, location, cook, degree), as.factor)) 

str(datos2)
'data.frame':   593 obs. of  24 variables:
 $ ID         : num  1 2 3 4 5 6 7 8 9 10 ...
 $ age        : num  20 45 20 18 18 18 18 18 18 18 ...
 $ gender     : Factor w/ 2 levels "0","1": 2 2 1 2 1 2 1 1 2 1 ...
 $ hiseilevel : Factor w/ 3 levels "1","2","3": 1 2 3 3 1 2 3 1 2 2 ...
 $ location   : Factor w/ 3 levels "1","2","3": 1 1 1 1 2 1 2 3 1 1 ...
 $ cook       : Factor w/ 2 levels "0","1": 1 2 1 1 2 1 1 2 1 1 ...
 $ degree     : Factor w/ 2 levels "0","1": 1 1 2 1 1 1 1 1 1 1 ...
 $ dairy      : num  3.42 6.64 3.42 1.5 1.26 ...
 $ olivesnuts : num  0.712 2 0.427 0.132 0 ...
 $ herbsspices: num  0.2795 2 0 0.0658 1.2137 ...
 $ fruit      : num  4.12 4.28 1.46 1.9 3.61 ...
 $ veget      : num  1.4 7.85 1.83 2.61 2.94 ...
 $ oliveoil   : num  0.4986 2.5014 2.5014 0.4986 0.0658 ...
 $ breadpasta : num  2.068 5.142 5.208 1.69 0.986 ...
 $ potato     : num  1.5 0 0.46 0.46 0.46 ...
 $ redmeat    : num  5.75 11.83 8.86 12.31 5.75 ...
 $ sweets     : num  6.367 1.956 18.871 8.899 0.921 ...
 $ whitemeat  : num  1.496 4.986 1.956 1.956 0.921 ...
 $ fish       : num  5.41 3.45 3.91 4.72 9.82 ...
 $ egg        : num  1.5 4.99 1.5 3.49 0.46 ...
 $ legumes    : num  3.452 5.408 2.877 2.301 0.921 ...
 $ alcoholicUA: num  0.21 0 0 0 3.38 ...
 $ fastfood   : num  0.329 0.263 0.921 0.46 0.707 ...
 $ precooked  : num  1.76 2.57 1.25 1.53 1.02 ...

Ahora, si queremos categorizar variables continuas en función de distintos niveles o cortes, podemos hacer lo siguiente:

# Utilizando condicionales if/else como en Excel:
datos <- openxlsx::read.xlsx("data/Dataset.xlsx", sheet=2)
datos$age_cat <- NA # Se inicializa la columna

for (i in 1:nrow(datos)) {
  if (datos$age[i] < 20) {
    datos$age_cat[i] <- "Muy joven"
  } else if (datos$age[i] >= 20 & datos$age[i] < 25) {
    datos$age_cat[i] <- "Menos joven"
  } else {
    datos$age_cat[i] <- "Ojo, que si haces deporte te lesionas"
  }
}

# Utilizando una función un poco más potente
datos <- openxlsx::read.xlsx("data/Dataset.xlsx", sheet=2)
datos$age_cat <- ifelse(datos$age < 20, 
                        "Muy joven", 
                        ifelse(datos$age >= 20 & datos$age < 25, 
                               "Menos joven", 
                             "Ojo, que si haces deporte te lesionas"))

# Utilizando funciones de dplyr

library(dplyr)
datos <- openxlsx::read.xlsx("data/Dataset.xlsx", sheet=2)
datos <- datos %>%
  mutate(age_cat = (if_else(
    age < 20,
    "Muy joven",
    if_else(age < 25, "Menos joven", "Ojo, que si haces deporte te lesionas")
  )))# Esto arroja un problema: solución
summary(datos$age_cat) # Es una variable de caracteres
   Length     Class      Mode 
      593 character character 
library(dplyr)
datos <- openxlsx::read.xlsx("data/Dataset.xlsx", sheet=2)
datos <- datos %>%
  mutate(age_cat = factor(
    if_else(
    age < 20,
    "Muy joven",
    if_else(
      age < 25,
      "Menos joven",
      "Ojo, que si haces deporte te lesionas")
    )
  )) 

# Otra forma
datos <- openxlsx::read.xlsx("data/Dataset.xlsx", sheet=2)
datos <- datos %>%
  mutate(age_cat = factor(case_when(
    age < 20 ~ "Muy joven",
    age >= 20 & age < 25 ~ "Menos joven",
    age >= 25 ~ "Ojo, te he dicho"
  )))

3. Tratamiento de bases de datos, análisis exploratorio y visualización: del tratamiento básico al {tidyverse} y a {ggplot2}

3.1 Tratamiento básico de bases de datos

Sin quererlo, antes ya hemos hecho algún tratamiento de datos. Incluso ya hemos utilizado parte del tidyverse. Ahora vamos a empezar a tratar los datos para conocerlos un poco mejor y usar algunas funciones básicas.

Ahora vamos a trabajar con un elemento (que creo) que es relevante y que nos permite trabajar con datos que estén en formato largo (long o narrow) o ancho (wide). En función del tipo de análisis que queramos hacer a veces tendremos que usar los datos en un formato o en otro. En el libro de The Epidemiologist R Handbook denominan a este proceso pivotar y recomiendan ver este capítulo de R for Data Science.

library(pacman)
pacman::p_load(rio, # Importar archivos
               here, # Localizar archivos
               kableExtra, # Crear y manipular tablas complejas
               tidyverse) # Paquete que incluye ggplot2, dplyr, tidyr, etc. Gestión de datos, como sabemos

# Importamos datos: https://github.com/appliedepi/epirhandbook_eng/raw/master/data/malaria_facility_count_data.rds 
library(DT)
count_data <- import("malaria_facility_count_data.rds")
datatable(count_data, options = list(pageLength = 5, filter = "top", lengthMenu = c(5, 10, 50, 100)))

¿Qué observamos? ¿Cómo son estos datos?

Estos datos podrían ser entradas de información desde distintas instalaciones sanitarias (columna newid como si fuese el número del informe). Es decir, cada fila es un “caso” (sujeto, observación para un ID) que incluye toda la información de este. A este tipo de datos se le considera en formato wide o ancho.

Por ejemplo, los datos que hemos trabajado antes sobre factores sociales y consumo de alimentos están en formato ancho.

Esta base de datos es, cito literal: "Each observation in this dataset refers to the malaria counts at one of 65 facilities on a given date, ranging from count_data$data_date %>% min() to count_data$data_date %>% max(). These facilities are located in one Province (North) and four Districts (Spring, Bolo, Dingo, and Barnard). The dataset provides the overall counts of malaria, as well as age-specific counts in each of three age groups - <4 years, 5-14 years, and 15 years and older."

De ancho a largo (wide-to-long)

# Pasamos de ancho a largo considerando, por ejemplo, la edad y el total

df_long <- count_data %>% 
  pivot_longer(
    cols = c(`malaria_rdt_0-4`, `malaria_rdt_5-14`, `malaria_rdt_15`, `malaria_tot`)
  )
# Forma 1 
df_long %>% 
  datatable(df_long, options = list(pageLength = 5, filter = "top", lengthMenu = c(5, 10, 50, 100)))
# Forma 2 con tydyselect helper
count_data %>% 
  pivot_longer(
    cols = starts_with("malaria_")
  )
# A tibble: 12,152 × 8
   location_name data_date  submitted_date Province District newid name    value
   <chr>         <date>     <date>         <chr>    <chr>    <int> <chr>   <int>
 1 Facility 1    2020-08-11 2020-08-12     North    Spring       1 malari…    11
 2 Facility 1    2020-08-11 2020-08-12     North    Spring       1 malari…    12
 3 Facility 1    2020-08-11 2020-08-12     North    Spring       1 malari…    23
 4 Facility 1    2020-08-11 2020-08-12     North    Spring       1 malari…    46
 5 Facility 2    2020-08-11 2020-08-12     North    Bolo         2 malari…    11
 6 Facility 2    2020-08-11 2020-08-12     North    Bolo         2 malari…    10
 7 Facility 2    2020-08-11 2020-08-12     North    Bolo         2 malari…     5
 8 Facility 2    2020-08-11 2020-08-12     North    Bolo         2 malari…    26
 9 Facility 3    2020-08-11 2020-08-12     North    Dingo        3 malari…     8
10 Facility 3    2020-08-11 2020-08-12     North    Dingo        3 malari…     5
# ℹ 12,142 more rows
# Forma 3 con posición de columnas (y hay otras formas)
count_data %>% 
  pivot_longer(
    cols = 6:9
  )
# A tibble: 12,152 × 8
   location_name data_date  submitted_date Province District newid name    value
   <chr>         <date>     <date>         <chr>    <chr>    <int> <chr>   <int>
 1 Facility 1    2020-08-11 2020-08-12     North    Spring       1 malari…    11
 2 Facility 1    2020-08-11 2020-08-12     North    Spring       1 malari…    12
 3 Facility 1    2020-08-11 2020-08-12     North    Spring       1 malari…    23
 4 Facility 1    2020-08-11 2020-08-12     North    Spring       1 malari…    46
 5 Facility 2    2020-08-11 2020-08-12     North    Bolo         2 malari…    11
 6 Facility 2    2020-08-11 2020-08-12     North    Bolo         2 malari…    10
 7 Facility 2    2020-08-11 2020-08-12     North    Bolo         2 malari…     5
 8 Facility 2    2020-08-11 2020-08-12     North    Bolo         2 malari…    26
 9 Facility 3    2020-08-11 2020-08-12     North    Dingo        3 malari…     8
10 Facility 3    2020-08-11 2020-08-12     North    Dingo        3 malari…     5
# ℹ 12,142 more rows

Nota que hemos pasado de 3.038 filas a 12.152 (ahora los datos son más largos). Concretamente, 4 veces más largos (3.038 *4 (3 grupo de edad + total), pero tienen menos columnas.

Este tipo de datos en formato largo se usan en ocasiones (depende de cada paquete) para hacer análisis de supervivencia. Los datos de panel también suelen estar en formato largo. En salud, podría ser el número de observaciones de un paciente en un hospital (ej: lo veo el lunes y veo cómo está, lo veo el martes y veo cómo está y cada fila representa un momento distinto)

De largo a ancho (long-to-wide)

Que el formato largo se relevante ocasiones, no indica que el formato ancho no lo sea. A veces, el formato ancho es más útil para hacer análisis. En el ejemplo de los factores sociales y demográficos y consumo de alimentos esa organización de datos nos permitió hacer modelos de regresión lineal múltiple.

Para transformar los datos de largo a ancho usaremos otros datos:

# Datos usados: https://github.com/appliedepi/epirhandbook_eng/raw/master/data/case_linelists/linelist_cleaned.rds

linelist <- import("linelist_cleaned.rds")
head(linelist)
  case_id generation date_infection date_onset date_hospitalisation
1  5fe599          4     2014-05-08 2014-05-13           2014-05-15
2  8689b7          4           <NA> 2014-05-13           2014-05-14
3  11f8ea          2           <NA> 2014-05-16           2014-05-18
4  b8812a          3     2014-05-04 2014-05-18           2014-05-20
5  893f25          3     2014-05-18 2014-05-21           2014-05-22
6  be99c8          3     2014-05-03 2014-05-22           2014-05-23
  date_outcome outcome gender age age_unit age_years age_cat age_cat5
1         <NA>    <NA>      m   2    years         2     0-4      0-4
2   2014-05-18 Recover      f   3    years         3     0-4      0-4
3   2014-05-30 Recover      m  56    years        56   50-69    55-59
4         <NA>    <NA>      f  18    years        18   15-19    15-19
5   2014-05-29 Recover      m   3    years         3     0-4      0-4
6   2014-05-24 Recover      f  16    years        16   15-19    15-19
                              hospital       lon      lat infector source wt_kg
1                                Other -13.21574 8.468973   f547d6  other    27
2                              Missing -13.21523 8.451719     <NA>   <NA>    25
3 St. Mark's Maternity Hospital (SMMH) -13.21291 8.464817     <NA>   <NA>    91
4                        Port Hospital -13.23637 8.475476   f90f5f  other    41
5                    Military Hospital -13.22286 8.460824   11f8ea  other    36
6                        Port Hospital -13.22263 8.461831   aec8ec  other    56
  ht_cm ct_blood fever chills cough aches vomit temp time_admission       bmi
1    48       22    no     no   yes    no   yes 36.8           <NA> 117.18750
2    59       22  <NA>   <NA>  <NA>  <NA>  <NA> 36.9          09:36  71.81844
3   238       21  <NA>   <NA>  <NA>  <NA>  <NA> 36.9          16:48  16.06525
4   135       23    no     no    no    no    no 36.8          11:22  22.49657
5    71       23    no     no   yes    no   yes 36.9          12:60  71.41440
6   116       21    no     no   yes    no   yes 37.6          14:13  41.61712
  days_onset_hosp
1               2
2               1
3               2
4               2
5               1
6               1
df_wide <- 
  linelist %>% 
  count(age_cat, gender)

df_wide
   age_cat gender   n
1      0-4      f 640
2      0-4      m 416
3      0-4   <NA>  39
4      5-9      f 641
5      5-9      m 412
6      5-9   <NA>  42
7    10-14      f 518
8    10-14      m 383
9    10-14   <NA>  40
10   15-19      f 359
11   15-19      m 364
12   15-19   <NA>  20
13   20-29      f 468
14   20-29      m 575
15   20-29   <NA>  30
16   30-49      f 179
17   30-49      m 557
18   30-49   <NA>  18
19   50-69      f   2
20   50-69      m  91
21   50-69   <NA>   2
22     70+      m   5
23     70+   <NA>   1
24    <NA>   <NA>  86
df_wide2 <- na.omit(df_wide) # Aprovecho y muestro cómo se eliminan na's

3.2 Análisis exploratorio de datos con {tidyverse} (among others) y {ggplot2}

Para este apartado se sigue la información de un capítulo del libro Fundamentos de ciencias de datos con R

Funcionamiento elemental de {ggplot2}

El paquete ggplot2nos permite crear gráficos que, honestamente, son una pasada. Una de las principales fuentes de información sobre el uso de este paquete lo podéis encontrar en el libro de Wickham y sus distintas ediciones (mira cuáles son las diferencias entre ellas, por ejemplo, la segunda incluía aspectos de Data Management).

De un curso reciente de Storytelling con R en la UCLM, me quedé con una imagen que presenta resumen clave de cómo funciona el paquete ggplot2 (que adaptaron de aquí, recurso también interesante):

Los gráficos se hacen de menos a más y se insertan distintas capas que ofrecen personalización. Hagamos un ejemplo sencillo también extraído de los materiales del mencionado curso para ver cómo funcionan las distintas capas:

library(ggplot2) # Recuerda que estos dos paquetes están en {tidyverse}
library(dplyr)  #cargamos los datos 'starwars'
library(ggthemes)
starwars <- starwars # redudante, pero para que se vea que se carga
# Paso 1: lienzo
ggplot(data = starwars)

# Paso 2: aesthethics mapping: eje x y eje y
ggplot(data = starwars,
       mapping=aes(x=height, y=mass))

# Paso 3: geometría
ggplot(data = starwars, 
       mapping=aes(x=height, y=mass))+
  geom_point() # Ojo outlier! ¿Quién es? Busca en la base de datos

# Paso 4: quitar outlier y representar de nuevo hasta aquí
starwars2 <- filter(starwars, name != "Jabba Desilijic Tiure")
ggplot(data = starwars2, 
       mapping = aes(x = height, y = mass)) +
  geom_point()

# Paso 5: distintas configuraciones
ggplot(data = starwars2,
       mapping = aes(x = height, y = mass, label = name)) +
  geom_point() +
  geom_text()

ggplot(data = starwars2,
       mapping = aes(x = height, y = mass, color = mass)) +
  geom_point(size = 2) +
  labs(
    title = "Star Wars: Altura y peso de los personajes",
    x = "Altura (cm)",
    y = "Peso (kg)",
    color = "Peso (kg)"
  )

ggplot(data = starwars2,
       mapping = aes(x = height, y = mass, color = gender)) +
  geom_point(size = 2) +
  labs(
    title = "Star Wars: Altura y peso de los personajes",
    subtitle = "por género",
    x = "Altura (cm)",
    y = "Peso (kg)",
    color = "Gender",
    caption = "Fuente: {dplyr}")+
  scale_colour_viridis_d()+ # Cambiamos la paleta de colores
  geom_abline(intercept=lm(mass~height, data=starwars2)$coefficients[1], slope=lm(mass~height, data=starwars2)$coefficients[2], color="darkblue")+  # Añadimos una línea de regresión
  theme_economist()

# Ojo: todo esto se puede hacer exactamente igual usando {dplyr}:

starwars2 %>% 
  ggplot(aes(x = height, y = mass, color=gender))+
  geom_point(size=2)+
  labs(
    title = "Star Wars: Altura y peso de los personajes",
    subtitle = "por género",
    x = "Altura (cm)",
    y = "Peso (kg)",
    color = "Gender",
    caption = "Fuente: {dplyr}")+
  scale_colour_viridis_d()+ # Cambiamos la paleta de colores
  geom_abline(intercept=lm(mass~height, data=starwars2)$coefficients[1], slope=lm(mass~height, data=starwars2)$coefficients[2], color="darkblue")+  # Añadimos una línea de regresión
  theme_economist()

En ese curso también nos mostraron un recurso web bastante útil: la galería de gráficos con R. Este recurso proporciona los MUCHÍSIMOS tipos de gráficos que se pueden crear en función de los datos que queremos analizar, así como sus definiciones y el código básico para obtenerlos (emulando una vignette o, mejor, una cheatsheet como esta). Pero sin duda, si no tenemos claro qué tipo de gráfico podemos obtener el recurso más relevante es el de from data to viz. Es más, el recurso anterior bebe de este último. Otro recurso similar puede encontrarse aquí.

Con todo, no podemos olvidar la siguiente máxima. Cuando queremos hacer un gráfico lo primero que tenemos que pensar es a qué audiciencia nos estamos dirigiendo ya que esto condiciona qué tipo de gráfico podemos usar: no es lo mismo una audiencia académica que la población general o que un decisor político; no es lo mismo dar una conferencia divulgativa que una científica. ¡Adáptate al contexto!

Let’s practice!

Para realizar el análisis exploratorio de datos elemental vamos a utilizar los datos sobre consumo de alimentos y factores sociales. Vamos a responder a estas tres cuestiones:

  1. ¿Cuál es la distribución del consumo de alimentos?
  2. ¿Cuál es el porcentaje hombres y mujeres que estudian o ciencias sociales o ciencias de la salud?
  3. ¿Cuál es la correlación entre el consumo de los distintos grupos de alimentos?

El exploratorio básico podría ser:

# Cargamos los datos

datos <- openxlsx::read.xlsx("data/Dataset.xlsx", sheet=2)

# Una función útil es la función summary() o explore_all()

summary(datos)
       ID           age            gender         hiseilevel       location    
 Min.   :  1   Min.   :17.00   Min.   :0.0000   Min.   :1.000   Min.   :1.000  
 1st Qu.:149   1st Qu.:18.00   1st Qu.:0.0000   1st Qu.:1.000   1st Qu.:1.000  
 Median :297   Median :19.00   Median :0.0000   Median :2.000   Median :1.000  
 Mean   :297   Mean   :20.21   Mean   :0.4199   Mean   :1.755   Mean   :1.583  
 3rd Qu.:445   3rd Qu.:21.00   3rd Qu.:1.0000   3rd Qu.:2.000   3rd Qu.:2.000  
 Max.   :593   Max.   :48.00   Max.   :1.0000   Max.   :3.000   Max.   :3.000  
      cook            degree           dairy          olivesnuts     
 Min.   :0.0000   Min.   :0.0000   Min.   : 0.000   Min.   :0.00000  
 1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.: 1.630   1st Qu.:0.06575  
 Median :0.0000   Median :0.0000   Median : 2.501   Median :0.21370  
 Mean   :0.2917   Mean   :0.2395   Mean   : 2.933   Mean   :0.34373  
 3rd Qu.:1.0000   3rd Qu.:0.0000   3rd Qu.: 3.707   3rd Qu.:0.42740  
 Max.   :1.0000   Max.   :1.0000   Max.   :14.507   Max.   :3.21370  
  herbsspices          fruit            veget            oliveoil     
 Min.   :0.00000   Min.   : 0.000   Min.   : 0.0000   Min.   :0.0000  
 1st Qu.:0.06575   1st Qu.: 1.321   1st Qu.: 0.9534   1st Qu.:0.4986  
 Median :0.27945   Median : 2.395   Median : 1.7068   Median :1.0000  
 Mean   :0.57098   Mean   : 2.895   Mean   : 2.2268   Mean   :1.1760  
 3rd Qu.:0.71233   3rd Qu.: 3.819   3rd Qu.: 2.7781   3rd Qu.:1.0000  
 Max.   :7.00000   Max.   :18.419   Max.   :26.9945   Max.   :6.0000  
   breadpasta        potato           redmeat           sweets      
 Min.   :0.000   Min.   : 0.0000   Min.   : 0.000   Min.   : 0.000  
 1st Qu.:1.282   1st Qu.: 0.4603   1st Qu.: 7.479   1st Qu.: 2.301  
 Median :1.918   Median : 1.4959   Median :10.912   Median : 4.871  
 Mean   :2.280   Mean   : 1.3234   Mean   :12.908   Mean   : 6.385  
 3rd Qu.:2.978   3rd Qu.: 1.4959   3rd Qu.:15.841   3rd Qu.: 7.978  
 Max.   :9.430   Max.   :17.5096   Max.   :80.989   Max.   :42.460  
   whitemeat           fish             egg            legumes      
 Min.   : 0.000   Min.   : 0.000   Min.   : 0.000   Min.   : 0.000  
 1st Qu.: 1.496   1st Qu.: 2.992   1st Qu.: 1.496   1st Qu.: 1.841  
 Median : 2.992   Median : 4.833   Median : 1.496   Median : 2.992  
 Mean   : 3.461   Mean   : 5.610   Mean   : 2.790   Mean   : 3.500  
 3rd Qu.: 3.951   3rd Qu.: 7.000   3rd Qu.: 3.490   3rd Qu.: 4.488  
 Max.   :33.005   Max.   :35.019   Max.   :31.510   Max.   :22.956  
  alcoholicUA         fastfood        precooked      
 Min.   : 0.0000   Min.   :0.0000   Min.   : 0.0000  
 1st Qu.: 0.0000   1st Qu.:0.3288   1st Qu.: 0.5260  
 Median : 0.2100   Median :0.5589   Median : 0.9041  
 Mean   : 0.7046   Mean   :0.6404   Mean   : 1.0286  
 3rd Qu.: 0.8400   3rd Qu.:0.9041   3rd Qu.: 1.3315  
 Max.   :19.1250   Max.   :2.8466   Max.   :12.2630  
datos %>% explore_all()

datos |> explore_all() # Lo mismo que arriba

explore_all(datos) # Lo mismo que arriba

# Otras estadísticas

mean(datos$legumes) # media de servings de legumbres
[1] 3.499822
sd(datos$fastfood) # desviación típica de servings de fastfood
[1] 0.4263786
# Como un summary más desarrollado y usando dplyr()

datos %>% select(dairy, fastfood, legumes) %>% descr()
Descriptive Statistics  
datos  
N: 593  

                     dairy   fastfood   legumes
----------------- -------- ---------- ---------
             Mean     2.93       0.64      3.50
          Std.Dev     2.02       0.43      2.75
              Min     0.00       0.00      0.00
               Q1     1.63       0.33      1.84
           Median     2.50       0.56      2.99
               Q3     3.71       0.90      4.49
              Max    14.51       2.85     22.96
              MAD     1.49       0.41      2.22
              IQR     2.08       0.58      2.65
               CV     0.69       0.67      0.79
         Skewness     1.84       1.23      2.71
      SE.Skewness     0.10       0.10      0.10
         Kurtosis     5.28       2.77     13.69
          N.Valid   593.00     593.00    593.00
        Pct.Valid   100.00     100.00    100.00

Para responder a la primera cuestión podemos usar histogramas u otro tipo de gráficos:

# 1. Distribución consumo alimentos: podríamos hacerlo de distintas formas e incluso mostrar varios histogramas a la vez en un mismo cuadro

ggplot(datos, aes(x = dairy)) +
  geom_histogram(binwidth = 0.5, fill = "#69b3a2", color = "white", alpha = 0.8) + # Personalización de colores
  labs(
    title = "Histograma de Lácteos",
    x = "Valores",
    y = "Frecuencia"
  ) +
  theme_minimal(base_size = 14) + # Tema minimalista
  theme(
    plot.title = element_text(hjust = 0.5, face = "bold", size = 16), # Centrar el título
    axis.title = element_text(face = "bold"), # Ejes con texto en negrita
    panel.grid.major = element_line(color = "gray85"), # Líneas de cuadrícula
    panel.grid.minor = element_blank() # Eliminar líneas menores
  )

# También podemos plantear gráficos de cajas 

ggplot(datos, aes(y = dairy)) +
  geom_boxplot(fill = "orange")

# Para trabajar con más de una variable hay que hacer los datos en formato largo, ¡ea!

datos_largos <- datos %>% 
  select(dairy, veget, legumes) %>% 
  pivot_longer(
    cols = c(dairy, veget, legumes)
  )
ggplot(datos_largos, aes(x=name, y = value)) +
  geom_boxplot(fill = "orange")

# Un poco feo y podríamos hacerlo más moderno y lustroso

ggplot(datos_largos, aes(x = name, y = value)) +
  geom_violin(color = "darkblue") +
  geom_boxplot(alpha = 0.5, 
               fill = "darkblue") +
  labs(title = "Consumo de alimentos diario",
       x = "Alimentos", 
       y = "Número de raciones medias") +
  theme_minimal()

Para estudiar el porcentaje de hombres y mujeres en los distintos tipos de estudios podemos estudiarlo a través de tablas de frecuencias elementales o usar algún gráfico de mosaico. Como siempre, se puede mejorar todo lo que queramos. Un ejemplo de esto:

# 2. Porcentaje de hombres/mujeres que estudian sociales/salud

table(datos$gender, datos$degree) # para tener una tabla sencilla de frecuencias
   
      0   1
  0 237 107
  1 214  35
prop.table(table(datos$gender,datos$degree),margin=1)
   
            0         1
  0 0.6889535 0.3110465
  1 0.8594378 0.1405622
# Los gráficos de mosaico son mu bonicos: uno básico
datos_fac = data.frame(gender=as.factor(datos$gender), degree=as.factor(datos$degree))
ggplot(data = datos_fac) +
  geom_mosaic(aes(x = product(degree, gender),
                  fill = degree))  

# Uno más completo:

# Instalar y cargar los paquetes necesarios
if (!require(ggplot2)) install.packages("ggplot2")
if (!require(ggmosaic)) install.packages("ggmosaic")

# Crear un conjunto de datos de ejemplo
library(ggmosaic)
set.seed(123)


# Gráfico de mosaico más "apañao"
ggplot(data = datos_fac) +
  geom_mosaic(aes(x = product(degree, gender), fill = degree), alpha = 0.9, color = "white") +
  scale_fill_viridis_d(option="D") +
  labs(
    title = "Distribución del tipo de estudios por género",
    subtitle = "Así más bonico",
    x = "Género",
    y = "Proporción",
    fill = "Tipo de estudios"
  ) +
  theme_minimal(base_size = 14) + # Tema limpio y profesional
  theme(
    plot.title = element_text(hjust = 0.5, size = 18, face = "bold"), # Título centrado
    plot.subtitle = element_text(hjust = 0.5, size = 14, face = "italic"), # Subtítulo
    axis.text.x = element_text(angle = 45, hjust = 1), # Rotar etiquetas del eje X
    legend.position = "right", # Leyenda a la derecha
    legend.title = element_text(face = "bold"), # Título de la leyenda en negrita
    panel.grid.major = element_blank(), # Sin cuadrícula principal
    panel.grid.minor = element_blank() # Sin cuadrícula secundaria
  )

  1. El estudio de la correlación de los distintos grupos de alimentos se puede abordar de varias formas:
# Correlación básica
cor(datos$fish, datos$redmeat, method="pearson")
[1] 0.2305928
cor.test(datos$fish, datos$redmeat, method="pearson") # Con significación

    Pearson's product-moment correlation

data:  datos$fish and datos$redmeat
t = 5.7611, df = 591, p-value = 1.346e-08
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 0.1529162 0.3054377
sample estimates:
      cor 
0.2305928 
# Matriz de correlaciones
datos_num <- datos[8:24]
mcor <- cor(datos_num) # matriz de correlaciones
round(cor(datos[8:24]),2) # con dos decimales
            dairy olivesnuts herbsspices fruit veget oliveoil breadpasta potato
dairy        1.00      -0.01       -0.04 -0.05  0.00    -0.01       0.08  -0.02
olivesnuts  -0.01       1.00        0.13  0.15  0.11     0.07       0.07   0.06
herbsspices -0.04       0.13        1.00  0.07  0.46     0.18       0.05   0.19
fruit       -0.05       0.15        0.07  1.00  0.20     0.12       0.01   0.02
veget        0.00       0.11        0.46  0.20  1.00     0.17       0.00   0.13
oliveoil    -0.01       0.07        0.18  0.12  0.17     1.00       0.13   0.05
breadpasta   0.08       0.07        0.05  0.01  0.00     0.13       1.00   0.04
potato      -0.02       0.06        0.19  0.02  0.13     0.05       0.04   1.00
redmeat      0.12       0.12        0.08 -0.04  0.05    -0.01       0.09   0.02
sweets       0.01       0.06        0.03  0.00 -0.06     0.02       0.07  -0.01
whitemeat    0.06      -0.01        0.10  0.08  0.10     0.03      -0.01   0.08
fish         0.07       0.15        0.16  0.19  0.15     0.06       0.05   0.08
egg          0.07       0.07        0.20  0.03  0.15     0.11       0.00   0.15
legumes     -0.05       0.14        0.12  0.21  0.18    -0.03       0.09   0.09
alcoholicUA  0.02       0.01        0.01  0.04 -0.05    -0.12      -0.06  -0.03
fastfood     0.00       0.04       -0.07 -0.12 -0.16    -0.11       0.09  -0.02
precooked    0.00       0.15        0.05  0.09  0.10    -0.01       0.02  -0.04
            redmeat sweets whitemeat  fish   egg legumes alcoholicUA fastfood
dairy          0.12   0.01      0.06  0.07  0.07   -0.05        0.02     0.00
olivesnuts     0.12   0.06     -0.01  0.15  0.07    0.14        0.01     0.04
herbsspices    0.08   0.03      0.10  0.16  0.20    0.12        0.01    -0.07
fruit         -0.04   0.00      0.08  0.19  0.03    0.21        0.04    -0.12
veget          0.05  -0.06      0.10  0.15  0.15    0.18       -0.05    -0.16
oliveoil      -0.01   0.02      0.03  0.06  0.11   -0.03       -0.12    -0.11
breadpasta     0.09   0.07     -0.01  0.05  0.00    0.09       -0.06     0.09
potato         0.02  -0.01      0.08  0.08  0.15    0.09       -0.03    -0.02
redmeat        1.00   0.05      0.10  0.23  0.05    0.05        0.05     0.21
sweets         0.05   1.00     -0.12 -0.07 -0.03   -0.03        0.04     0.19
whitemeat      0.10  -0.12      1.00  0.23  0.21    0.08        0.04    -0.10
fish           0.23  -0.07      0.23  1.00  0.05    0.06        0.03    -0.05
egg            0.05  -0.03      0.21  0.05  1.00    0.05        0.03    -0.02
legumes        0.05  -0.03      0.08  0.06  0.05    1.00       -0.04     0.00
alcoholicUA    0.05   0.04      0.04  0.03  0.03   -0.04        1.00     0.09
fastfood       0.21   0.19     -0.10 -0.05 -0.02    0.00        0.09     1.00
precooked      0.05   0.08      0.05  0.12  0.02    0.05        0.01     0.17
            precooked
dairy            0.00
olivesnuts       0.15
herbsspices      0.05
fruit            0.09
veget            0.10
oliveoil        -0.01
breadpasta       0.02
potato          -0.04
redmeat          0.05
sweets           0.08
whitemeat        0.05
fish             0.12
egg              0.02
legumes          0.05
alcoholicUA      0.01
fastfood         0.17
precooked        1.00
# Gráficos de correlación que pueden ser útiles para explorar los datos 
corrplot(mcor, method="circle", diag=F) # se puede arreglar tanto como queramos

# Este también es mu bonico, pero aquí tenemos muchos datos y habría que arreglarlo estéticamente
datos_num %>% 
  ggpairs(progress = FALSE)+
  theme(plot.margin = unit(c(1, 1, 1, 1), "cm"))

Para rematar un gráfico de barras que no tiene nada que ver:

# Gráfico de barras (en este caso recuento de malaria)

ggplot(data = df_long) +
  geom_col(
    mapping = aes(x = data_date, y = value, fill = name),
    width = 1
  )

4. Introducción a los análisis estadísticos aplicados a la Economía y al Análisis de Políticas Públicas

Para saber qué tipo de análisis debemos aplicar, es importante conocer los tipos de datos (tipo de variable) tenemos y cuáles se adaptan mejor a cada modelo. Con el tiempo y la experiencia, las transformaciones de datos son comunes. Esto no significa que tengamos que torturar los datos hasta que obtengamos el resultado que intuimos podría ser adecuado: al contrario. Antes de plantear cualquier tipo de modelo estadístico tenemos que tener muy claro el diseño de nuestra investigación y qué han hecho investigaciones con temáticas y/o diseños similares. Esto nos permitirá saber el estado de la cuestión, lo común y lo novedoso que podemos aportar. Las transformaciones de datos o la creación de variables proxies surgirán por distintas razones como, por ejemplo, alguna limitación en el diseño. (casi) Siempre que trabajemos con honestidad y transparencia, obtendremos resultados satisfactorios.

Con este cuadro resumen (preparado por una antigua y buena profesora que me dio estadística cuando cursaba el máster e imagino que inspirado en el libro de Bioestadística amigable, porque el libro lo conocí en el máster) podemos ver distintos tipos de análisis en función de nuestros datos:

Recuerda también que hay paquetes que te permiten importar a R directamente las bases de datos desde las webs y navegar desde el programa. Esto es tremendamente útil. Un ejemplo sería el paquete {eurostat}. Otro ejemplo, usado en finanzas, sería descargar datos desde Yahoo! Finance.

# Eurostat

library(eurostat)
eurostat::eu_countries
# A tibble: 27 × 3
   code  name     label                                           
   <chr> <chr>    <chr>                                           
 1 BE    Belgium  Belgium                                         
 2 BG    Bulgaria Bulgaria                                        
 3 CZ    Czechia  Czechia                                         
 4 DK    Denmark  Denmark                                         
 5 DE    Germany  Germany (until 1990 former territory of the FRG)
 6 EE    Estonia  Estonia                                         
 7 IE    Ireland  Ireland                                         
 8 EL    Greece   Greece                                          
 9 ES    Spain    Spain                                           
10 FR    France   France                                          
# ℹ 17 more rows
# Yahoo! Finance

# Descarga de datos -------------------------------------------------------

inicio <- as.Date("2000-01-01") #fecha de inicio
final <- as.Date(Sys.Date())  #fecha de final
activos <- c("^IBEX") #activos que quiero cargar
getSymbols(activos,
source = "yahoo",
from = inicio,
to = final) #extrae los datos
[1] "IBEX"
datos_ibex <- IBEX[, 2:4]
# Gráfico con ggplot2
library(ggplot2)
datos_ibex2 <- data.frame(fecha = index(datos_ibex), coredata(datos_ibex))
ggplot(datos_ibex2, aes(x = fecha)) +
geom_line(aes(y = IBEX.High, color = "Cierre"), size = 1) +
geom_line(aes(y = IBEX.Low, color = "Bajo"), size = 1) +
geom_line(aes(y = IBEX.Close, color = "Alto"), size = 1) +
labs(
title = "Evolución del IBEX35",
x = "Fecha",
y = "Precio",
color = "Leyenda"
) +
scale_color_manual(values = c(
"Cierre" = "darkblue",
"Bajo" = "green3",
"Alto" = "darkred"
)) +
theme_minimal() +
theme(
legend.position = "top",
legend.title = element_text(face = "bold"),
axis.text.x = element_text(angle = 45, hjust = 1)
)

# Gráfico dinámico --------------------------------------------------------

# Esto toma las columnas de máx, mín y cierre
dygraph(datos_ibex,
xlab = "Periodo",
ylab = "Cotización",
main = "Cotización IBEX 35 del 2000 a la actualidad") %>% dyRangeSelector()

Con todo, para trabajar algunos modelos estadísticos, usaremos varias bases de datos ya utilizadas.

4.0. Diferencias en medias, proporciones, etc.

Estos análisis son importantes y suelen incluirse en las tablas de descriptivos de cualquier trabajo académico. Nos dan información relevante y muchas veces ya nos hablan de su relevancia como factores asociados (estudios transversales) o relacionados (en estudios longitudinales) en modelos un poco más complejos.

t.test(datos$redmeat~datos$gender, var.equal=T) # t-Student. Si no son normales U Mann-Whitney

    Two Sample t-test

data:  datos$redmeat by datos$gender
t = -2.7066, df = 591, p-value = 0.006994
alternative hypothesis: true difference in means between group 0 and group 1 is not equal to 0
95 percent confidence interval:
 -3.5164750 -0.5591071
sample estimates:
mean in group 0 mean in group 1 
       12.05248        14.09027 
var.test(datos$redmeat~datos$gender) # Si son o no varianzas iguales

    F test to compare two variances

data:  datos$redmeat by datos$gender
F = 1.0018, num df = 343, denom df = 248, p-value = 0.9926
alternative hypothesis: true ratio of variances is not equal to 1
95 percent confidence interval:
 0.7930661 1.2600626
sample estimates:
ratio of variances 
           1.00184 
chisq.test(datos$gender, datos$degree) # Para estudiar diferencias en frecuencias en tablas de contingencia

    Pearson's Chi-squared test with Yates' continuity correction

data:  datos$gender and datos$degree
X-squared = 22.126, df = 1, p-value = 2.554e-06

4.1 Regresión lineal y logística simple

Recuerda que antes de hacer cualquier regresión hay que estudiar cuáles son las condiciones necesarias para poder hacer incluir las variables en una regresión. En general, las regresiones se utilizan para predecir (mejor R^2), pero también es común usarlas para estudiar asociaciones (sin fines de predicción), donde el poder de explicación de la varianza no tiene por qué ser tan alto.

En la regresión lineal sus supuestos básicos son: - Linealidad - Independencia del error o residuos: no hay correlación entre errores - Normalidad de error: error normalmente distribuido - Homocedasticidad (ausencia de heterocedasticidad ;-) ): hace referencia a que la varianza de los errores de la regresión es la misma para cada observación (así la variabilidad de la respuesta será la misma para valores bajo y altos de la variable independiente)

# Lineal

model <- lm(datos$redmeat ~ datos$veget)
summary(model) # model es un objeto con info dentro! Abre con $

Call:
lm(formula = datos$redmeat ~ datos$veget)

Residuals:
    Min      1Q  Median      3Q     Max 
-14.552  -5.487  -1.793   2.823  67.569 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  12.5250     0.5100  24.557   <2e-16 ***
datos$veget   0.1721     0.1560   1.103     0.27    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 9.095 on 591 degrees of freedom
Multiple R-squared:  0.002055,  Adjusted R-squared:  0.0003661 
F-statistic: 1.217 on 1 and 591 DF,  p-value: 0.2704
# Asunciones del modelo
plot(model, which=c(1:6))

par(mfrow=c(2,3))
# Linealidad en residuals vs fitted, no hay patrones
# Independencia del error 
library(lmtest)
dwtest(model) # Aprox a 2 no hay autocorrelación - Durbin Watson

    Durbin-Watson test

data:  model
DW = 1.9958, p-value = 0.4785
alternative hypothesis: true autocorrelation is greater than 0
# Normalidad del error: el gráfico Q-Q. PUntos cerca de la línea.

# Homocedasticidad: varianza de residuos constante. Con gráfico si resid vs fitted no tiene patrón - Si, por ejemplo forma de embudo, heterocedasticidad. El test: Breusch-Pagan o White
bptest(model) # BP - Sale heterocedasticidad - Se puede corregir aplicando transformaciones logarítmicas en las variables

    studentized Breusch-Pagan test

data:  model
BP = 4.361, df = 1, p-value = 0.03677
bptest(model, ~fitted(model)+I(fitted(model)^2)) # White: no asume relacion lineal entre residuos y variables independ.

    studentized Breusch-Pagan test

data:  model
BP = 4.3984, df = 2, p-value = 0.1109
model <- lm(log(datos$redmeat ~ log(datos$veget)))
Error in log(datos$redmeat ~ log(datos$veget)): Argumento no numérico para una función matemática
# Logística: esta la vemos por encimísima:
model2 <- glm(datos$degree~datos$sexo, family="binomial") # hay más familias, depende de los datos
Error in model.frame.default(formula = datos$degree ~ datos$sexo, drop.unused.levels = TRUE): tipo (NULL) inválido para la variable 'datos$sexo'

Da error, ¿por qué? Soluciónalo con lo que ya sabemos.

datos$degree_f <- as.factor(datos$degree)
datos$gender_f <- as.factor(datos$gender)
model2 <- glm(datos$degree_f~datos$gender_f, family="binomial")
model2 <- glm(degree_f~gender, data=datos, family="binomial")
summary(model2)

Call:
glm(formula = degree_f ~ gender, family = "binomial", data = datos)

Coefficients:
            Estimate Std. Error z value Pr(>|z|)    
(Intercept)  -0.7952     0.1165  -6.828 8.62e-12 ***
gender       -1.0154     0.2164  -4.693 2.69e-06 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 652.84  on 592  degrees of freedom
Residual deviance: 628.70  on 591  degrees of freedom
AIC: 632.7

Number of Fisher Scoring iterations: 4
OR <- exp(coef(model2)) # OR y sus IC
IC <- exp(confint(model2))
data.frame(OR,IC)
                   OR    X2.5..   X97.5..
(Intercept) 0.4514768 0.3580582 0.5655336
gender      0.3622587 0.2344588 0.5486563

4.2 Regresión lineal múltiple

Para realizar la regresión lineal múltiple se utiliza información y bases de datos de un manual titulado Introductory Econometrics: A modern approach y las soluciones están aquí. El ejemplo sirve para preguntarse sobre el efecto de distintas variables en la cantidad de manzanas ecológicas demandadas por una familia, incluyéndose entre las variables independientes el precio de las manzanas ecológicas, el de las normales y otras de caracter socioeconómico. En este ejemplo se cumple la teoría básica de la demanda: efecto negativo del precio de la manzana eco y efecto positivo del precio de la manzana normal.

library(wooldridge)
apple <- wooldridge::apple

model_app <- lm(data=apple, ecolbs~ ecoprc+regprc+faminc+hhsize+educ+age)
summary(model_app)

Call:
lm(formula = ecolbs ~ ecoprc + regprc + faminc + hhsize + educ + 
    age, data = apple)

Residuals:
   Min     1Q Median     3Q    Max 
-2.570 -1.146 -0.595  0.515 39.940 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  1.056770   0.892650   1.184    0.237    
ecoprc      -2.861237   0.591991  -4.833 1.68e-06 ***
regprc       3.006077   0.712308   4.220 2.79e-05 ***
faminc       0.002195   0.002865   0.766    0.444    
hhsize       0.063094   0.067780   0.931    0.352    
educ         0.034313   0.045314   0.757    0.449    
age          0.001389   0.006763   0.205    0.837    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 2.486 on 653 degrees of freedom
Multiple R-squared:  0.04022,   Adjusted R-squared:  0.0314 
F-statistic: 4.561 on 6 and 653 DF,  p-value: 0.0001519
# Sacado también de https://herioscarlanda.wordpress.com/wp-content/uploads/2018/10/wooldridge-2009-introduccic3b3n-a-la-econometrc3ada-un-enfoque-moderno.pdf

4.3 Diferencias en diferencias (diff-in-diff)

El análisis diff-in-diff o DiD es útil para experimentos naturales o quasiexperimentales (randomizados es complicado). Este tipo de análisis nos permite estudiar diferencias entre un grupo de intervención y otro de control en función del tiempo, cuestiones geográficas o por valores tras la adopción de una política o tras la sucesión de algún suceso relevante. Permite comparación transversal (entre grupo de control e intervención) y antes-después (cómo estaba cada grupo antes-después del suceso). Es un método popular para hacer inferencia causal en microeconometría aplicada.

4.3.0. Introducción

Para estudiar este método usaré el desarrollo del modelo del manual Learning Microeconometrics with R de Adams (2021). Las medidas repetidas en un individuo permite infererir causalmente el efecto del tratamiento y descubrir características no observadas, siendo esto útil en conjuntos de datos de panel.

Considera que los datos son como los de la figura que aparece más abajo.

Estamos interesadxs en el efecto de X en Y que se denota por b. Este efecto varía por cada individuo (faltaría un sub “i”, vaya). Considera que observamos dos Xs diferentes para el mismo individuo: dos niveles de educación distinto. Así, podemos estimar el efecto de la educación del individuo en el resultado de interés (Y). Podemos medir la diferencia en la renta asociado a un aumento en la educación: el efecto causal, en otras palabras, de la educación sobre la renta. Sin embargo, si observamos el mismo individuo con dos niveles de educación distintos, estamos observando al individuo en dos veces diferentes. Esas caractéristicas no observadas (U) también pueden cambiar entre el tiempo: ha un efecto de la educación y del tiempo en Y. El efecto del tiempo se denota por d. Esta imagen muestra que la relación entre X e Y es confusa. Hay una backdoor relationship entre X e Y, que funciona a través del tiempo: puedes ir de X a Y o de X a tiempo, de tiempo a U y de U a Y. Si el efecto de Tiempo a X y a U es 1, la regresión de X a Y es b+d. Para estimar b necesitamos una forma de estimar d.

4.3.1 Primeras diferencias

Observamos un outcome de interés para el individuo i en dos momentos del tiempo, 1 y 2: \[ y_{it} = a + b x_{it} + d v_{it} + \epsilon_{it} \] donde el outcome está determinado por la variable política/intervención “x” y el efecto no observado del tiempo por “u”. Tomando primeras diferencias podemos considerar el efecto no observado del tiempo (u):

\[ y_{i2} - y_{i1} = b (x_{i2} - x_{i1}) + d (v_{i2} - v_{i1}) + \epsilon_{i2} - \epsilon_{i1} \]

Si asumimos que el efecto del tiempo es 1, se puede añadir un intercepto para medir el efecto del tiempo. Lo desarrollamos en el siguiente ejemplo

Ejemplo con datos de panel simulados y primera estimación

Considera estos datos de panel simulados:

set.seed(123456789)
N <- 1000
T <- 2
a <- -2
b <- 3
d <- -4
u <- c(5,6) #características no observadas
x1 <- runif(N)
x2 <- x1+runif(N)
x <- as.matrix(rbind(x1,x2)) # matriz TxN
e <- matrix(rnorm(N*T), nrow=T) # otras características no observadas adicionales
y <- a + b*x + d*u + e

Estos datos permiten observar el cambio de x para un individuo. Dado que estamos interesados en el efecto de este cambio en x podemos mirar primeras diferencias. Hacemos regresión del cambio en X sobre el cambio en Y.

diffy <- y[2, ] - y[1, ]
diffx <- x[2, ] - x[1, ]
lm1 <- lm(diffy ~ diffx)
summary(lm1)

Call:
lm(formula = diffy ~ diffx)

Residuals:
    Min      1Q  Median      3Q     Max 
-5.1334 -0.9315  0.0577  0.8596  4.7124 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -3.77099    0.09013  -41.84   <2e-16 ***
diffx        2.74810    0.15416   17.83   <2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.415 on 998 degrees of freedom
Multiple R-squared:  0.2415,    Adjusted R-squared:  0.2408 
F-statistic: 317.8 on 1 and 998 DF,  p-value: < 2.2e-16

Si se usa el intercepto, b se aproxima a 3 y el intercepto (d) a -4.

4.3.2 Diff-in-diff

Ahora vamos a desarrollar un ejemplo para aplicar un diff-in-diff siguiendo los pasos proporcionados en este tutorial (el ejemplo también está desarrollado en el libro que estamos siguiendo).

Los datos del ejemplo de Card and Krueger (1994) consisten en unas encuestas en restaurantes antes y después de que el estado de New Jersey aumentase el salario mínimo. El objetivo era estudiar el impacto en el empleo. En 1991, era de $4.25/h y en 1992, de $5.05/h. Así, se podría comparar simplemente el impacto en el empleo de la variación del salario mínimo tras la aplicación normativa (antes-después). Sin embargo, no parece suficiente el antes-después sin más porque otros factores podrían haber intervenido. Por tanto, se utiliza como referencia el estado vecino de Pensilvania como grupo de control (o sin tratamiento), considerando que que los restaurantes de ese estado no sufrieron impacto por el cambio normativo en New Jersey. La teoría económica ortodoxa indica que subidas de salario mínimo conllevan pérdidas de empleo, ¿será así?

Empezamos explorando los datos ¿sabemos ya cómo hacerlo?:

library(openxlsx)
rest <- read.xlsx("data/datos_diff.xlsx")
viridis_color <- viridis(5)[1] # saco un color
summary(rest)
    bonus              chain             co_owned              date      
 Length:820         Length:820         Length:820         Min.   :33553  
 Class :character   Class :character   Class :character   1st Qu.:33919  
 Mode  :character   Mode  :character   Mode  :character   Median :33921  
                                                          Mean   :33924  
                                                          3rd Qu.:33925  
                                                          Max.   :33969  
                                                          NA's   :410    
     empft            emppt          firstinc         hrsopen     
 Min.   : 0.000   Min.   : 0.00   Min.   :0.0000   Min.   : 7.00  
 1st Qu.: 2.000   1st Qu.:11.00   1st Qu.:0.1500   1st Qu.:12.00  
 Median : 6.000   Median :17.00   Median :0.2000   Median :15.00  
 Mean   : 8.239   Mean   :18.75   Mean   :0.2214   Mean   :14.45  
 3rd Qu.:12.000   3rd Qu.:25.00   3rd Qu.:0.2500   3rd Qu.:16.00  
 Max.   :60.000   Max.   :60.00   Max.   :1.0000   Max.   :24.00  
 NA's   :18       NA's   :14      NA's   :123      NA's   :11     
    inctime         meals               ncalls          nmgrs       
 Min.   : 2.00   Length:820         Min.   :0.000   Min.   : 0.000  
 1st Qu.:13.00   Class :character   1st Qu.:0.000   1st Qu.: 3.000  
 Median :19.00   Mode  :character   Median :1.000   Median : 3.000  
 Mean   :20.08                      Mean   :1.447   Mean   : 3.452  
 3rd Qu.:26.00                      3rd Qu.:2.000   3rd Qu.: 4.000  
 Max.   :52.00                      Max.   :9.000   Max.   :10.000  
 NA's   :97                         NA's   :249     NA's   :12      
     nregs          nregs11      observation             open       
 Min.   :1.000   Min.   :0.000   Length:820         Min.   : 0.000  
 1st Qu.:3.000   1st Qu.:2.000   Class :character   1st Qu.: 6.500  
 Median :3.000   Median :3.000   Mode  :character   Median : 7.000  
 Mean   :3.601   Mean   :2.689                      Mean   : 8.076  
 3rd Qu.:4.000   3rd Qu.:3.000                      3rd Qu.:10.500  
 Max.   :8.000   Max.   :7.000                      Max.   :11.500  
 NA's   :28      NA's   :39                         NA's   :11      
     pctaff          pentree           pfry            psoda      
 Min.   :  0.00   Min.   :0.410   Min.   :0.6700   Min.   :0.410  
 1st Qu.: 15.00   1st Qu.:0.940   1st Qu.:0.8500   1st Qu.:1.000  
 Median : 50.00   Median :1.030   Median :0.9400   Median :1.050  
 Mean   : 48.87   Mean   :1.338   Mean   :0.9315   Mean   :1.045  
 3rd Qu.: 80.00   3rd Qu.:1.843   3rd Qu.:1.0100   3rd Qu.:1.090  
 Max.   :100.00   Max.   :3.950   Max.   :1.3700   Max.   :1.490  
 NA's   :454      NA's   :36      NA's   :45       NA's   :30     
    region              sheet         special             state          
 Length:820         Min.   :  1.0   Length:820         Length:820        
 Class :character   1st Qu.:119.0   Class :character   Class :character  
 Mode  :character   Median :237.5   Mode  :character   Mode  :character  
                    Mean   :246.5                                        
                    3rd Qu.:372.0                                        
                    Max.   :522.0                                        
                                                                         
    status              type              wage_st          emptot     
 Length:820         Length:820         Min.   :4.250   Min.   : 0.00  
 Class :character   Class :character   1st Qu.:4.500   1st Qu.:14.50  
 Mode  :character   Mode  :character   Median :5.000   Median :20.00  
                                       Mean   :4.806   Mean   :21.03  
                                       3rd Qu.:5.050   3rd Qu.:25.50  
                                       Max.   :6.250   Max.   :85.00  
                                       NA's   :41      NA's   :26     
    pct_fte     
 Min.   : 0.00  
 1st Qu.:15.38  
 Median :32.65  
 Mean   :34.03  
 3rd Qu.:54.05  
 Max.   :90.36  
 NA's   :32     
rest %>% 
  explore_all()

rest <- rest %>% filter(!is.na(rest$wage_st))
tapply(rest$wage_st, list(rest$state, rest$observation), mean)
             February 1992 November 1992
New Jersey        4.612134      5.080849
Pennsylvania      4.630132      4.617465
rest_before <- subset(rest, observation=="February 1992") # Saco la información de febrero de 1992 y muestro la función subset()
rest_after <- subset(rest, observation=="November 1992") # Saco la información de noviembre de 1992

rest_before %>% ggplot(aes(x=wage_st))+ # Histograma de salario
  geom_histogram(fill=viridis_color, color="black", bins=20)+
  labs(title="Salario en ambos estados Feb 1992", x="Salario", y="Frecuencia")+
  xlim(3,7)+
  #ylim(0,350)+
  geom_vline(xintercept=c(4.25,5.05), color="black")

rest_after %>% ggplot(aes(x=wage_st))+ # Histograma de salario
  geom_histogram(fill=viridis_color, color="black", bins=20)+
  labs(title="Salario en ambos estados Nov 1992", x="Salario", y="Frecuencia")+
  xlim(3,7)+
 # ylim(0,350)+
  geom_vline(xintercept=c(4.25,5.05), color="black")

Una vez explorados los datos, procedemos a estudiar las diferencias de medias ode empleo entre grupos (primeras diferencias:

differences <- rest %>%
  group_by(observation, state) %>%
  summarise(emptot = mean(emptot, na.rm = TRUE))
differences
# A tibble: 4 × 3
# Groups:   observation [2]
  observation   state        emptot
  <chr>         <chr>         <dbl>
1 February 1992 New Jersey     20.4
2 February 1992 Pennsylvania   23.4
3 November 1992 New Jersey     21.5
4 November 1992 Pennsylvania   21.6
# Treatment group (NJ) before treatment
njfeb <- differences[1,3]

# Control group (PA) before treatment
pafeb <- differences[2,3]

# Treatment group (NJ) after treatment
njnov <- differences[3,3]

# Control group (PA) after treatment
panov <- differences[4,3]

Ahora sacamosel Efecto Medio del Tratamiento (Average Treatment Effect, ATT):

(njnov-njfeb)-(panov-pafeb) # Diferencia entre la diferencia de noviembre y febrero en cada estado
    emptot
1 2.876178
(njnov-panov)-(njfeb-pafeb) # Diferencia entre la diferencia de los estados en noviembre y febrero
    emptot
1 2.876178

Por último, podemos obtener el estimador DiD:

# Arreglamos los datos

rest <- mutate(rest, time=ifelse(observation=="November 1992", 1,0),
               treated=ifelse(state=="New Jersey", 1,0))
did_model <- lm(emptot ~ time + treated + time:treated, data = rest)
summary(did_model) # Nota que los resultados no coinciden con el tutorial porque eliminamos NAs antes

Call:
lm(formula = emptot ~ time + treated + time:treated, data = rest)

Residuals:
    Min      1Q  Median      3Q     Max 
-15.929  -6.419  -0.972   4.528  64.581 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)    23.429      1.079  21.714   <2e-16 ***
time           -1.823      1.542  -1.183   0.2374    
treated        -3.010      1.203  -2.503   0.0125 *  
time:treated    2.876      1.714   1.678   0.0938 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 9.282 on 756 degrees of freedom
  (19 observations deleted due to missingness)
Multiple R-squared:  0.00892,   Adjusted R-squared:  0.004987 
F-statistic: 2.268 on 3 and 756 DF,  p-value: 0.07932

El coeficiente de la itneracción tiempo-tratamiento es el estimador DiD. El tratamiento (aumento de salario mínimo) lleva en media a un incrempento de empleo en Nueva Jersey de 2.88 FTE (que es el número de trabajadores a tiempo completo (incluyendo gestores) más 0.5 veces el número de trabajadores a tiempo parcial)

4.4 Otros ejemplos de técnicas estadísticas aplicadas a economía y políticas públicas: análisis de supervivencia y transiciones entre estados

Estos ejemplos os los enseño desde mi propio PC pues son propios y los datos no son públicos.

5. Reproducibilidad con Quarto

Quarto es un sistema de publicación técnico y científico de código abierto. Permite crear documentos reproducibles e interactivos para distintos lenguajes: Python, R, Julia u Observable JS. En otras palabras, multilenguaje. Es la siguiente generación de R Markdown en Posit (RStudio) con otras características adicionales. Permite exportar en distintos formatos, pero con htmlya hay bastante juego.

Es el programa con el que he desarrollado esta hoja de ruta y os mostraré lo elemental para que podáis crear contenido similar. Lo (que creo que es más) importante:

  1. Configurar el YAML: basta con copiarlo de otro más desarrollado y modificarlo, probando qué pasa. Lo básico es título, fecha y autoría.
  2. Las almohadillas (#) son importantes para dividir e indexar el contenido.
  3. Los bloques de código se abren con tres acentos graves y el lenguaje de programación y se cierran con otros tres acentos graves (estos acentos son el opuesto al acento/tilde en castellano).
  4. Configurar bien la carpeta donde tenemos el proyecto pues se convertirá en el directorio de trabajo y podremos acceder fácilmente. Evitar rutas excesivamente largas porque si no puede fallar al leer las rutas de imágenes o bases de datos.
  5. Es otro lenguaje que funciona bastante bien con IA generativa. Haz las preguntas adecuadas y testea.
  6. Si te animas, únelo con github y comparte lo que desarrolles :).

Declaración de uso de IA generativa (no sé si hace falta ponerlo, pero creo que es lo correcto)

Para realizar este documento se han utilizado ChatGPT y Copilot (integración con GitHub) para autogenerar partes de código, para detección de errores, para corregirlos y para limpiar el código.